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19. ABSTRACT (continued) 


“^In contrast with time domain amplitudes, the elastic structure of the travel path causes 
little variability in spectral amplitudes. When the mantle Q is constrained to a value 
of approximately 100 at depths greater than 120 km, the inversion for Q and moment pro¬ 
duces moments that remain constant with distance. Based on the best models available, 
surface waves from NTS explosions should be larger than surface waves from East Kazakh 
exp osiops with the same moment. Estimated scalar moments for the largest East Kazakh 
explosions since 1976 are smaller than the estimated moments for the largest NTS explo- 
sions for the same time period. This contrasts dramatically with m. for the same events 
which show the opposite behavior. ^ b 
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I. INTRODUCTION 


Accurate and unbiased surface wave measurements are required if 
surface waves are to be used for treaty monitoring purposes. 
Conventional time domain M s measurements provide a rough estimate of 
the size of the seismic source. Differences in travel path, however, 
cause variations in M g that are unacceptable in cases where 0.1 
magnitude unit may be important. On short paths (less than about 30 
degrees), the seismogram is dominated by the Airy phase and 20 second 
amplitudes are difficult or impossible to measure accurately. On long 
paths, the surface waves are welI-dispersed but attenuation 
differences can cause large amplitude variations. Since explosion 
yield estimation is usually done by defining a magnitude versus yield 
curve at one or more test sites, path differences and differences in 
station coverage can lead to a biased result when using this curve to 
estimate yields at an uncalibrated test site. 

In this report, we describe a method for obtaining more stable 
and reliable surface wave measurements by using all of the information 
contained in the surface wave to estimate the Green’s function for the 
path. We use Rayleigh waves from explosions at the Nevada Test Site 
(NTS) and the Soviet East Kazakh test site to find Green’s functions 
for the paths between these test sites and SRO and WWSSN stations. 
The phase of the Green’s function can be used to construct a phase- 
matched filter, which may be used to recover spectral amplitudes for 
surface waves from other events observed along nearby paths. The 
moment of the explosion can then be found by dividing the spectral 
amplitude by the amplitude of the Green’s function for the path. 

Figure 1 illustrates the difficulty in measuring surface wave 
amplitudes using standard methods. The top figure is a Rayleigh wave 
from an NTS explosion recorded at Station GOL, at a distance of 977 
km. At this close distance many frequencies arrive at the same time, 
so that the seismogram is dominated by the Airy phase. A reliable 20 
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Seismogram and spectrum, NTS-GOL 


This figure illustrates the difficulty in measuring surface wave 
amplitudes using standard methods. The top figure is a Rayleigh wave 
from the NTS explosion SANDREEF recorded at Station GOL. The 
seismogram is dominated by the Airy phase, so a reliable 20 second 
amplitude cannot be measured. The lower figure shows the average 
spectrum (solid line) and * 1 standard deviation (dotted line) for 24 
spectra at Station GOL. The amplitude depends strongly on the 
frequency measured. 











second amplitude cannot be measured. Measuring the spectral amplitude 
instead of the seismogram removes this problem. The lower figure 
shows the average spectrum (solid line) and ♦_ 1 standard deviation 
(dotted line) for 24 spectra at Station GOL. Another problem still 
remains, however. The amplitude depends strongly on the frequency 
measured, changing by a factor of 2 between 0.035 Hz and 0.05 Hz, and 
by a factor of 5 between 0.02 Hz and 0.05 Hz. Figure 2 shows the 
amplitude and phase of the spectrum after the path correction has been 
applied. Notice that in addition to correcting for path effects, the 
path correction is also a frequency correction, flattening the 
spectrum over the frequency band of interest. 

The procedures recommended in this report are based on our 
experience applying surface wave analysis methods to real and 
synthetic data over the past few years. Bache, et al. (1978a) used 
inversion of surface wave dispersion to estimate the path structure 
between NTS and WWSSN Stations ALQ and TUC. Bache, et al. (1978b) 
used these structures to estimate the source size of NTS explosions. 
Rodi, et al. (1978) inverted Rayleigh waves from explosions at the 
French Sahara test site to estimate path structure in Africa. Our 
procedure is similar to that used by Bache, et al. (1978a), except 
that we use phase-matched filtering to recover phase velocities and 
spectral amplitudes; we invert for the Q structure of the path as well 
as the elastic structure; and the programs have been stream)ined so 
many paths can be processed in an efficient manner. Tests of the 
phase-matched filtering and inversion programs are given by Wang, et 
at. (1981). The computer programs used in this study are described in 
detail in a report by Stevens, et al. (1982). 
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Path corrected amplitude and phase for the path NTS-GOL. The variation in amplitude with 
frequency is greatly reduced. The phase correction can be used to determine phase shifts 
(such as Rayleigh wave reversals) between events. 











II. RAYLEIGH WAVES FROM AN EXPLOSION SOURCE 


Our objective is to find the Green’s function for a given source 
to receiver path which can be used to determine the explosion moment. 
Our approach to the problem has been to assume that the average path 
structure can be adequately modeled by a plane-layered, attenuating 
elastic medium. Since the explosions are located in a common source 
region, but the Rayleigh waves travel along different paths, we also 
use an approximation which allows us to separate the spectrum into 
source region and path dependent functions (Bache, et aL., 1978a). 
This approximation assumes energy conservation and no mode conversion 
when a Rayleigh wave crosses a boundary between adjacent pIane-Iayered 
structures. The accuracy of this approximation has been demonstrated 
for selected cases by comparison with finite difference calculations 
by Glover and Harkrider (1981). We also assume an isotropic explosion 
source, neglecting non isotropic source effects such as explosion- 
induced block motion and tectonic strain release. In the construction 
of the path Green’s functions, we try to use events with little 
evidence of non isotropic source effects. 

Given these assumptions, we can write the vertical component of 
the Rayleigh wave as a product of functions each of which depends only 
on the source, the source region, or the path: 

u(w,r) * M Q 'S x («)S 2 («)exp(~ 7 2 r)/Ja # sin(r/a e ) «xp[i(* 0 -wr/c 2 )](1) 


where Mq' is the normalized moment of the explosion discussed below, 
r is the source to receiver distance, a e is the radius of the earth, 
C 2 is the frequency dependent phase velocity for the path structure, 
72 is the frequency dependent attenuation coefficient due to the 
anelastic structure of the path, and Iq is the initial phase which is 
equal to -3r/4 for an explosion with step function moment. Sj(«) and 
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S 2 (w) are real functions that depend on the earth structure of the 
source region and path respectively. Sj(w) is positive for a shallow 
source. Using the notation of Harkrider (1964, 1970), 



E} and E 3 are the eigenfunctions of horizontal displacement and normal 
traction respectively, and are written in Harkrider’s notation as Ej = 
Us(h)/£g, ^3 = a s(b)/(&o/ c ) 0 * * s * hear modulus at the source). 
They are functions of the source depth h and the large scale 
properties of the earth structure, cj and C 2 are the phase velocities 
of the source and path structures, and Apj and Ap 2 are the frequency 
dependent Rayleigh wave amplitude functions for the source and path 
regions which depend on the large scale properties of the earth 
structure. 

In the following analysis, the source region elastic structure 
is assumed to be constrained by independent geophysical evidence and 
the same structure will be used for all events at a given test site. 
The quantities which must be determined in Equation ( 1 ) are C 2 , S 2 , 
and 72 . 

The normalized moment Mg is equal to 3/J^/a^Mg, where Mg is the 
conventional definition of the moment and a and p are the 
compressionaI velocity and shear velocity at the source. For a 
Poisson solid Mq a Mg. We define Mg this way so that it will be a 
measure of the strength of surface wave excitation by the source, and 
so the Green’s function will be independent of local material 
properties. With this definition, two nearby explosions in media with 
different Poisson’s ratios will produce surface waves of the same 
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amplitude if they have the same Mq, provided that the change in local 
material properties is not large enough to significantly change the 
long period eigenfunctions in Equation 2. 
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III. EXPLOSION DATA FOR THE NEVADA TEST SITE AND EAST KAZAKH 


The data used in this study were Rayleigh wave seismograms from 
underground nuclear explosions at the Nevada Test Site (NTS) and the 
Soviet East Kazakh test site (EKZ). Stations were selected which had 
good signal-to-noise ratios for at least two events. The location of 
the recording stations relative to the test sites are shown in Figures 
3, 4, and 5. Twenty-two explosions at the East Kazakh test site were 
recorded at 12 SRO stations. Forty-eight NTS explosions were recorded 
at 12 SRO stations and 24 WWSSN stations. 

The WWSSN data were hand digitized by ENSCO Incorporated, and 
installed at the Center for Seismic Studies in Arlington, Virginia. 
Except for occasional timing errors, we found the WWSSN data to be of 
very good quality, and we have had little difficulty in applying the 
inversion procedure to this data. All of the data processing was done 
on the VAX 11/780 at the Center for Seismic Studies. 

The paths shown in Figures 3 through 5 cover a variety of earth 
structures and a considerable range of distances. The paths from NTS 
to WWSSN stations are relatively short, ranging from 900 to 4400 km 
and are all continental paths. The SRO stations are all quite distant 
from NTS, ranging from 5600 to 14,000 km. Three of these paths are 
almost entirely oceanic, while the remainder are mixed continental and 
oceanic. Most of the paths from East Kazakh to SRO stations are 
between 2000 and 5000 km, but two of the paths are longer than 10,000 
km. Seismograms are nearly identical for all events observed along 
the same path, except for anomalous events with reversed Rayleigh 


waves. 




























IV. DATA ANALYSIS AND INVERSION PROCEDURE 


A block diagram of the method used for the data processing is 
shown in Figure 6. We start by finding a set of Rayleigh wave 
seismograms from a given source region and recorded at the same 
station which show consistent signals and have good signaI-to-noise 
ratios. We then find the phase velocities, group velocities and 
amplitude spectra for each of these seismograms using narrow band 
filtering (Archambeau, et aL ., 1966), Dziewonski, et aL., 1969, 
Savino, et at., 1980) and phase-matched filtering (Herrin and Goforth, 
1977; Wang, et aL., 1981). Group velocities are found first by 
applying a set of narrow band filters at successive center frequencies 
and computing the arrival times of the energy envelopes. The phase 
velocities are found next by integrating the group delays to find the 
phase as a function of frequency, forming a phase-matched filter, 
cross-correlating the filter with the seismograms, windowing the 
compressed signal, and unwrapping the residual phase spectrum to 
correct the initial phase estimate. We have found that this procedure 
produces very accurate and consistent measurements of phase velocity. 
The residual phase is then differentiated to obtain an improved group 
velocity estimate, and the spectral amplitudes of the Rayleigh waves 
are found by taking the Fourier transform of the windowed, compressed 
signaI. 

After applying this procedure to several seismograms, we check 
the results to make sure there are no phase-reversed seismograms or 
other inconsistencies in the dispersion curves. We then average the 
phase and group velocities and spectral amplitudes for the selected 
seismograms. This averaging serves three functions - it smooths the 
data, reduces the chance of error from inaccurate times or locations, 
and provides standard deviations on the data which are used as weights 
by the inversion code. 
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Figure 6. Block diagram of surface wave analysis procedure. 





















Quality control during path correction processing is very 
important, and is the most time consuming part of the procedure. The 
following procedure was used to select the data used for path 
correction processing. First, the data for each station were plotted 
and low signal-to-noise ratio data were eliminated (Figure 7). The 
selected data where replotted and marked at the 3 km/sec group 
velocity point on the seismograms to help identify timing errors in 
the data (Figure 8). This step was essential for the hand-digitized 
WWSSN data. Ten percent of the data was found to contain timing 
errors. Most of these are minute errors, however, and for these 
seismograms we simply adjusted the start time by 60 seconds. More 
troublesome were smaller errors of 10 or 20 seconds. Since we could 
not easily guess the error in these seismograms, they were not used in 
the data processing. 

Since it was necessary to process a large amount of data («• 500 
seismograms for ~ 50 paths), the interactive procedure described by 
Stevens, et aL. (1982) was partially automated. A trial set of input 
parameters for the data analysis code TELVEL were selected. These 
parameters are listed in Table 1. The variable parameters are the 
frequency band and filter Q for narrow-band filtering, the windowing 
interval for phase-matched filtering, output frequencies, and 2ir 
branch for phase velocity. The phase velocity branch was found by 
estimating the phase velocity at one frequency and finding the branch 
closest to this value. The values listed in Table 1 were based on our 
experience with digital data (where the correct branch can often be 
found by the long-period behavior), on worldwide average values from 
Oliver (1962) and on phase velocities measured in the western United 
States by McEvilly (1964). One seismogram for each path was processed 
interactively to make sure the input parameters were appropriate. 

After the initial processing it is necessary to select the most 
consistent results and average them for inversion. To check this, we 
plot all of the phase and group velocity curves to identify 
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inconsistent data. There is almost always some bad data resulting 
from interference, digitization errors, phase reversals, incorrect 
phase velocity branches, or timing errors. In the upper figure in 
Figure 9, we show phase and group velocities for all seismograms 
processed at Station LON. Three seismograms were rejected on the 
basis of the figure. The phase velocity curve which is offset from 
the rest resjlted from a 10 second timing error. The two curves which 
are inconsistent at long-periods are the result of digitization 
errors. This is a common problem in hand-digitized data when there is 
a strong Airy phase - narrow-band filtering finds the Airy phase even 
at the longest periods, so a nearly constant group velocity curve 
results. We have never seen this effect in digital data. The lower 
figure shows the average phase and group velocities together with 
_♦ one standard deviation curves for the remaining nine seismograms at 
Station LON. 

In Figure 10, we show the final phase and group velocities found 
for NTS to WWSSN stations, NTS to SRO stations, and East Kazakh to SRO 
stations. The phase and group velocities are very consistent for 
different events observed along a given path. Typical phase velocity 
standard deviations are 0.1 to 0.2 percent, while group velocity 
standard deviations are a few tenths of a percent. The standard 
deviations of all quantities increase near the ends of the observable 
frequency band, and also near spectral dips. 
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Figure 9. 


Phase and group velocities for path from NTS to LON (1155 km), top 
figure shows velocities for each seismogram. Bottom figure shows 
average (solid) and * 1 a curves after deleting data with timing and 
digitization errors. Origin id numbers are marked for inconsistent 
events. 
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gure 10. Phase and group velocities for paths from NTS to’WWSSN stations 
(top), NTS to SRO stations (middle) and East Kazakh to SRO stations 
(bottom) Solid lines are phase velocities; dashed are group 

velocities. The three paths with high phase velocity in the NTS to 
SRO figure are predominantly oceanic paths. The three dispersion 
curves with very low group velocities at 0.03 Hz in the East Kazakh 
figure correspond to the paths with thick, low velocity crusts to the 
south of the test site. 
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V. INVERSION FOR SHEAR VELOCITY STRUCTURE 


The shear velocity structure is found by simultaneously 
inverting the phase and group velocities. The method used is similar 
to the method used by Bache, et at. (1978a). The inversion program 
attempts to fit the data with a smooth model by minimizing the 
integral of ld/7/dzl over the structure, except across discontinuities 
at depths specified by the user. The inversion method is iterative, 
with the inversion model for each iteration becoming the starting 
model for the next iteration. The final results are completely 
independent of the starting model, and the initial starting model is 
necessary only to allow the calculation of partial derivatives of 
phase and group velocities with respect to model parameters. 

Inversion is performed for shear velocity only, as there is not 
enough information in the fundamental mode phase and group velocities 
to determine the compressionaI velocities and densities at the same 
time. We use data on crust and upper mantle properties from Oobrin 
(1976) to constrain these parameters. Compressions! velocity is 
constrained to be consistent with a Poisson’s ratio of 0.27, while 
density is constrained by a Birch’s law given by p = 0.65 p ♦ 400 (in 
MKS units). 

A 20 layered model is used for the inversion. The thickness of 
the layers and depth of the model are derived from the frequency 
content of the data. The maximum depth is approximately equal to the 
longest wavelength c/f, and layers are thickest near the maximum depth 
and become thinner near the free surface. The inversion is damped, so 
that the model chosen has approximately five or six degrees of freedom 
depending on the quality of the data. This allows an excellent fit to 
the data in most cases, and produces a smooth model free of the 
spurious oscillations which would result from allowing more degrees of 
f reedom. 











Discontinuities in the structure are explicitly allowed at the 
crust-mantle boundary and in some cases at shallower depths to allow 
for near surface sedimentary layers. The depths of the 
discontinuities are found in the following way. We first perform an 
inversion to find the best model with no discontinuities. Then we 
insert discontinuities at each layer boundary, invert, and search for 
minima in the resulting data misfit. Discontinuities are inserted at 
these depths, and a final inversion is performed with the 
discontinuities incIuded. 

In Figure 11, we show the data misfit minima found using the 
initial structure (with no discontinuities) shown in Figure 12. Data 
misfits are shown for several degrees of freedom. The data misfit is 
normalized to the standard deviations on the data, so a data misfit of 
1.0 means that the average misfit between the calculated and observed 
phase and group velocities is one standard deviation. Minima are 
found close to 7.5 km and 40 km, so discontinuities are inserted at 
these depths. The final model is also shown in Figure 12. 
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VI. INVERSION FOR R STRUCTURE AND SCALAR MOMENT 


The final step in the inversion procedure is the inversion of 
the average spectral amplitude for R and moment. The average spectral 
amplitude is the linear average of the spectra of the best events 
recorded for a given path, defined by 


“M ■ S 


N 

ZZu.oo 

i=i 


( 3 ) 


where u; are the individual spectra and N is the number of 
frequencies. The relative moments of the individual spectra may be 
found using the least squares relation 


= IZu («)«:(«)/ 2Zu(w)u(tt) (4) 

u u 

so once the absolute moment corresponding to the average spectrum is 
determined, the moments of all the events are determined. Standard 
deviations on the normalized spectra (spectra divided by their 
relative moments) are used as weights in the inversion. The standard 
deviations are defined by 

ff 2 (u(w)) = j-y TZ |«00 (5) 


In Figure 13, we show the spectra of individual seismograms and the 
average spectrum with standard deviations for the path from NTS to 
Station LON. 
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Rewriting Equation ( 1 ), we define the data used for the 
inversion to be the ratio between the synthetic spectral amplitude 
S (u) (for unit moment, without including attenuation) and the observed 
amplitude spectrum u(u): 


*'"!©■ - r '• “o' 


S(w) is defined by 


S(w) = Sj(«)S 2 (e»)/Ja^ sin '(r/.J 


The attenuation coefficients 7 (w) can be expanded in terms of the Q 
structure of the path using the phase velocity partial derivatives 
(Mitchell, 1975; Anderson, et at., 1965). Q is assumed to be depth 
dependent, but frequency independent. Explicitly, assuming no bulk 
losses 



4 9c 
3 a ; Oa. Q. 


( 8 ) 


where the index i indicates the layer in the inversion model and N is 
the total number of layers. a,, p\, and Q-, are the compressionaI 
velocity, shear velocity and Q, respectively, in each layer. 

In the Q/moment inversion, the amplitude of the spectral ratio 
in Equation ( 6 ) determines the moment, while the shape determines the 
Q structure. We use a smoothness criterion, as in the inversion for 
shear velocity structure to damp oscillations in the solution. Since 
there is some correlation between low velocity zones and low Q zones 
in the earth, we invert for p/ Q, and find the smoothest /9/Q structure 




consistent with the data. The Q/moment inversion is less robust than 
inversion for structure, since there is an additional free parameter 
(the moment) and since spectra are less smooth than dispersion curves. 
As a result, the solution must be more heavily damped by allowing 
fewer degrees of freedom in the inversion. The minimum number of 
degrees of freedom is two - one to determine the moment and one to 
find the average /J/Q value for the path. With exactly two degrees of 
freedom, Q in each layer is proportional to the shear velocity, with 
the proportionality constant determined by the inversion. More 
degrees of freedom allow more shape in the Q structure. 

We impose one additional constraint on the Q inversion. Tests 
on synthetic seismograms show that while the Q inversion procedure 
accurately reproduces the Q structure of the crust, lower Q zones 
beneath the crust are not well resolved, and the Q at depth is usually 
overestimated. This causes the moments measured at distant stations 
to be underestimated. We can reduce this problem by including an 
estimated value for the average mantle Q at depth as a data point. 
With this constraint, the Q inversion is stable and physicaIly 
reasonable Q structures are obtained with an inversion that allows 
approximately 3.5 degrees of freedom in the model. 







VII. SHEAR VELOCITY INVERSION RESULTS 


Several representative shear velocity structures obtained by 
inversion are shown in Figure 14. Three continental and one 
predominantly oceanic structures are shown for NTS. The continental 
structures are all quite similar. The path to SHA has unusually low 
velocities in the upper few kilometers. Note the corresponding low 
group velocities for this path in Figure 10. The path from East 
Kazakh to KONO crosses a stable, relatively fast shield structure. 
The path to MAJO is similar, but somewhat slower and more complex. 
The path to SHIO crosses the Tibetan plateau, and the inferred 
structure has a thick, low velocity crust. The path to ANTO has low 
velocities in the upper part of the crust. 

Our main interest here is in the effect of these structures on 
surface wave amplitudes In Figure 15, we show the path amplitude 
functions S 2 calculated using all of the path structures. The elastic 
structure of the source-to-receiver path has surprisingly little 
effect on the spectral amplitudes. The maximum variation in log (S 2 ) 
is only 0.15 at 20 seconds. The paths with the largest values of S 2 
are the thick crust, low velocity paths, while the smallest values 
occur along oceanic paths. 

These results mean that the surface wave amplitudes are not 
changed very much by the elastic properties of the travel path. This 
is, of course, only true for spectral amplitudes, not for time domain 
amplitudes. Differences in dispersion can cause large differences in 
time domain surface wave amplitudes, especially at distances less than 
30 degrees. On long or complex paths, when a plane-layered model is 
not a good approximation, interference will cause more spectral 
amplitude variation than is predicted by the smooth functions in 
Figure 15. 




















Path Amplitude Function ($2)- All Paths 



Figure 15. Path amplitude function Sg for all paths from NTS and East Kazakh. 











VIII. q CONSTRAINTS AND INVERSION RESULTS 


To invert for the path Q structures and the explosion moments, 
we need the excitation function Sj for the source region structures, 
and an estimate of mantle Q at depth to be used as a constraint. 

The Soviet test site is located in the Kazakh Fold System near 
latitude 50 degrees north and longitude 79 degrees east. According to 
Beliayevsky, et aL. (1968), the results of a deep seismic sounding 
survey show that the average crustal thickness in the area is 50 km. 
Close to the test site, however, the crustal thickness decreases to 
about 45 km. Compressional velocities in the crust are very similar 
to velocities in the stable shield and platform regions of the 
Eurasian continent, except for the absence of low velocity surface 
sediments in the Kazakh Fold System. Grand, et aL. (1984) show that 
the shear velocity in the crust near the source region is relatively 
slow, similar to the shear velocity of the platforms to the south and 
east, but slower than the velocity along the shield paths to Europe. 

For our model of the East Kazakh source region, we use a 
modification of the structure inferred for the path to MAIO. This 
path crosses the Turanian Plateau southwest of the test site. The 
crustal thickness of this structure is 45 km. We have increased the 
inferred velocities in the top 15 km (see Figure 16) to remove the 
influence of the low velocity sediments in the plateau. The structure 
in the upper crust is then consistent with velocities from 
Beliayevsky, et aL. (1968). We have reduced the velocity in the top 
1.5 km to 5.0 km/sec (compressional) and 2.7 km/sec (shear) with a 
density of 2.7 gm/cm^. These values are appropriate for weathered 
hard rock such as granite. 
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Figure 16. Source region structure used for the Esst Kazakh test site. 

The structure found for East Kazakh to MAIO (dashed line) 
was modified to construct the East Kazakh structure. 








The models used for NTS and East Kazakh are listed in Tables 2 
and 3 and shown in Figure 17. The NTS model is based on the NTS-TUC 
crustal structure of Bache, et at. (1978a), merged with the upper 
mantle structure of Anderson and Hart (1976, 1978). The shallow 
structure (upper 4 km) is the Pahute Mesa structure from Bache, et aL. 
(1975). The shear velocity at the source depth of 1 km is 2000 m/sec, 
and the shear velocity in the top 500 m is 1000 m/sec. 

The excitation functions for the two source regions are shown in 
Figure 18. The excitation function for NTS is greater than the 
excitation function for East Kazakh at all frequencies of interest, 
and particularly at higher frequencies. We therefore expect an 
explosion at NTS to generate larger surface waves than an explosion at 
East Kazakh with the same moment. 

We also need an estimate of worldwide average mantle Q to use as 
a constraint on the inversions, so we have searched through the 
literature to find a representati ve set of Q estimates. Anderson and 
Hart (1978) find a worldwide average Q of 90 to 110 between 80 and 300 
km depth. Canas and Mitchell (1978) find a minimum mantle Q of 100 to 
200 in the Pacific Ocean with the minimum at about 150 km depth, while 
Canas and Mitchell (1981) find a minimum mantle Q of 60 to 150 in the 
Atlantic Ocean at 100 to 150 km depth. Cara (1981) finds a mantle Q 
of 80 to 100 between 100 and 300 km depth for paths across the Pacific 
Ocean. Lee and Solomon (1979) find very low mantle Q values of about 
20 in the Western United States between 60 and 160 km, 25 in the east- 
central United States below 135 km, and 70 in the Central Pacific 
below 60 km. Singh (1982) finds a mantle Q of about 50 in southern 
Asia. Patton (1980) finds mantle Q values of 60 to 150 in western 
Eurasia at depths below 50 km, 100 to 150 in eastern Eurasia at depths 
below 100 km, and 300 to 500 in the stable platforms of Eurasia at 
depths below 50 km. There is a considerable amount of uncertainty in 
all of these results, and considerable variation in the estimates. 
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TABLE 2 

NTS SOURCE REGION STRUCTURE 


Thickness P-Velocity S-Velocity Density 

(km) (km/sec) (km/sec) (9/cm 3 ) 


0.5 2.00 1.00 1.70 

1.0 3.30 2.00 2.10 

1.5 4.50 2.70 2.40 

1.0 5.90 3.40 2.75 

8.0 5.96 3.52 2.78 

9.0 6.11 3.61 2.80 

10.0 6.37 3.76 2.84 

14.0 7.90 4.42 3.20 

20.0 8.05 4.50 3.30 

15.0 8.10 4.50 3.30 

.0 8.00 4.40 3.30 

.0 7.90 4.30 3.25 

.0 7.90 4.30 3.25 

40.0 8.00 4.40 3.30 

• 8.50 4.70 3.50 











TABLE 3 


EAST KAZAKH SOURCE REGION STRUCTURE 


Thickness 

(km) 

P-Velocity 

(km/sec) 

S-Velocity 

(km/sec) 

Density 
( 9 /cm 5 ) 

1.5 

5.00 

2.70 

2.70 

1.0 

5.40 

3.00 

2.70 

3.0 

5.90 

3.30 

2.70 

11.0 

6.10 

3.40 

2.70 

5.5 

6.31 

3.54 

2.70 

5.5 

6.60 

3.70 

2.81 

5.5 

6.74 

3.78 

2.86 

5.6 

6.78 

3.81 

2.88 

6.5 

6.80 

3.81 

2.88 

8.0 

8.15 

4.57 

3.37 

9.4 

8.14 

4.57 

3.37 

10.9 

8.11 

4.55 

3.36 

12.8 

8.06 

4.53 

3.34 

15.0 

8.05 

4.52 

3.34 

17.5 

8.07 

4.53 

3.35 

20.4 

8.12 

4.56 

3.36 

23.9 

8.15 

4.58 

3.38 

27.9 

8.16 

4.58 

3.38 

32.6 

8.15 

4.58 

3.37 

38.2 

8.12 

4.56 

3.36 

m 

8.10 

4.55 

3.36 









Source Region Elaetlc Structure 
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In the inversions described in this paper, we set /J/Q equal to 
40 with a standard deviation of 1 between 120 and 150 km depth. Since 
the shear velocity at this depth is between 4000 and 5000 m/sec, Q is 
constrained to be slightly more than 100 in this depth range. The 
inversion algorithm treats this "constraint" as a data point and tries 
to find a Q structure that both fits the observed spectral amplitude 
data and has the preferred Q value in the specified depth range. The 
Q value from the inversion may vary from the preferred value if the 
data requires it. 

Several representative Q structures are shown in Figure 19. 
Inferred crustal Q is quite variable, with very low values for some 
complex paths such as East Kazakh to ANTO in Turkey. Seismograms 
recorded at KONO (from East Kazakh) are rich in high frequencies, 
which may imply a high crustal Q along this path. 

In Figure 20, we show the inferred attenuation coefficients for 
all paths. In Table 4, we have listed the average 20 second values of 
attenuation coefficients, together with the structural parameters, for 
paths from NTS to WWSSN stations and SRO stations, for the three 
oceanic paths out of NTS, for paths from East Kazakh to SRO stations 
and for all paths. Differences in attenuation may be responsible for 
large amplitude variations. Rayleigh waves observed along paths with 
attenuation coefficients of 1.0 x 1.0" 7 /meters and 3.0 x 10~ 7 /meters 
(the * one standard deviation values for all paths) will differ by 
0.09 magnitude units at 1000 km, 0.43 magnitude units at 5000 km, and 
0.87 magnitude units at 10,000 km. 

















Attenuation Coefficients (y)-AII Paths 



Attenuation coefficients for all paths froai NTS and East Kazakh. 











TABLE 4 

AVERAGE 20 SECOND VALUES FOR NTS AND 
EAST KAZAKH PATHS (MKS UNITS) 





IX. ESTIMATED MOMENTS FOR NTS AND EAST KAZAKH EXPLOSIONS 


We have estimated explosion moments for the 40 NTS explosions 
and 16 East Kazakh explosions which were recorded at two or more 
stations. The estimated moments (in MKS units, newton-meters), 
together with their network standard deviations and the number of 
stations used in the average a^e given in Tables 5 and 6. The yields 
given in Table 5 are from Marshall, et a L. (1979). The network 
standard deviations in (log) moment are quite small, about 0.1. Even 
for recent NTS explosions, which include data from several distant SRO 
stations, the standard deviations are only about 0.15. 

As a check on the results, we calculated residual station 
corrections for all paths. These are listed in Tables 7 and 8. These 
station corrections are designed to minimize the variance in moment 
over the network at stations observing each event. Most of the 
station corrections are small, particularly for shorter paths. 
Eighteen of the 24 WWSSN stations have station residuals of 0.1 or 
less. The station corrections increase with distance from the source, 
which is to be expected because of the increasing complexity of the 
travel paths, but the largest corrections for some very complex paths 
are only about 0.3. 

Sykes and Cifuentes (1984) measured M s values for a similar set 
of East Kazakh explosions. Our moment estimates are very consistent 
with their reported M s values. Using the ten normal events that the 
two data sets have in common we find the following moment/M s relation: 

log M 0 = M s ♦ 11.86 * 0.03 

where the * 0.03 is the standard deviation over all ten events. 

It is interesting to compare the results for recent events with 
theoretical explosion models. Since April 1976, underground 
explosions have been limited to a yield of 150 kilotons by the 
Threshold Test Ban Treaty. In Figure 21, we show the estimated 
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TABLE 5 


NTS EXPLOSIONS 

Number 


Event 

Date 

Area* 

of 

Stations 

Log 

Moment 

Standard 

Deviation 

Remarks 

Corduroy 

12/03/65 

Y 

5 

16.10 

0.14 


Dumont 

05/19/66 

Y 

5 

15.98 

0.13 


Piledriver 

06/02/66 

C 

6 

15.84 

0.12 

Yield 62 KT, Reversed 

Commodore 

05/20/67 

Y 

5 

16.27 

0.26 

Yield 250 KT 

Scotch 

05/23/67 

P 

7 

16.46 

0.18 

Yield 155 KT 

Zaza 

09/27/67 

Y 

6 

16.32 

0.14 


Lanpher 

10/18/67 

Y 

4 

15.87 

0.06 


Stinger 

03/22/68 

P 

3 

16.17 

0.09 


Sled 

08/29/68 

P 

10 

16.22 

0.09 


Noggin 

09/06/68 

Y 

6 

16.05 

0.06 


Purse 

05/07/69 

P 

9 

16.42 

0.06 


Ti jeras 

10/14/70 

Y 

5 

16.01 

0.13 


Oscuro 

09/21/72 

Y 

3 

15.99 

0.11 


Almendro 

06/06/73 

P 

10 

16.78 

0.11 


Escabosa 

07/10/74 

Y 

7 

16.14 

0.15 


Portmanteau 

08/30/74 

Y 

2 

16.08 

0.00 


Topgallant 

02/28/75 

Y 

2 

15.83 

0.09 


Tybo 

05/14/75 

P 

10 

16.50 

0.10 


Sti(ton 

06/03/75 

P 

5 

16.17 

0.12 


Mi zzen 

06/03/75 

Y 

5 

16.02 

0.15 


Mast 

06/19/75 

P 

12 

16.58 

0.15 


Cheshire 

02/14/76 

P 

4 

16.79 

0.39 


Estuary 

03/09/76 

P 

4 

16.79 

0.05 


Pool 

03/17/76 

P 

3 

16.52 

0.09 


Strait 

03/17/76 

Y 

3 

16.28 

0.15 


Marsi1ly 

04/05/77 

Y 

2 

16.01 

0.06 

150 KT Limit Begins 

Sandreef 

11/09/77 

Y 

15 

16.27 

0.15 


Farallones 

12/14/77 

Y 

2 

16.22 

0.02 


Iceburg 

03/23/78 

Y 

2 

16.15 

0.07 


Lowba11 

07/12/78 

Y 

4 

15.96 

0.07 


Pan i r 

08/31/78 

P 

8 

15.97 

0.15 


Rummy 

09/27/78 

Y 

16 

16.12 

0.15 


Farm 

12/16/78 

P 

7 

15.89 

0.18 


Pepato 

06/11/79 

P 

17 

16.07 

0.13 


Hearts 

09/06/79 

Y 

16 

16.12 

0.15 


Sheepshead 

09/26/79 

P 

14 

15.95 

0.14 


Pyramid 

04/16/80 

Y 

5 

15.65 

0.08 


Colwick 

04/26/80 

P 

7 

15.91 

0.18 


Kash 

06/12/80 

P 

8 

16.16 

0.09 


Taf i 

07/25/80 

P 

7 

16.10 

0.09 



*Y = Yucca Flat, P = Pahute Mesa, C = Climax Stock 
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TABLE 6 


EAST KAZAKH EXPLOSIONS 


Number 


Date 

of 

Stations 

Log 

Moment 

Standard 

Deviation 

Remarks 

07/04/76 

3 

15.73 

0.18 


06/11/78 

5 

15.94 

0.04 


09/15/78 

6 

15.74 

0.10 


11/29/78 

9 

15.75 

0.11 


06/23/79 

10 

15.83 

0.09 


08/04/79 

10 

15.94 

0.12 


08/18/79 

9 

15.69 

0.20 

Reversed 

10/28/79 

6 

15.85 

0.08 


12/02/79 

9 

15.96 

0.11 


12/23/79 

3 

15.66 

0.20 

Large body waves 

09/14/80 

10 

15.86 

0.09 

Reversed 

10/12/80 

7 

15.94 

0.09 


04/22/81 

7 

15.92 

0.09 


09/13/81 

4 

15.87 

0.32 

Very low at KAAO 

10/18/81 

3 

15.91 

0.13 


12/27/81 

5 

15.92 

0.09 


04/25/82 

5 

15.85 

0.10 


12/05/82 

5 

15.92 

0.12 







TABLE 7 


NTS SURFACE WAVE PATHS 


Number Residual 

Average of Statior 

Distance (km) Events Correctic 


















TABLE 8 

EAST KAZAKH SURFACE WAVE PATHS 


Station 

Average 
Distance (km) 

Number 

of 

Events 

Residua1 
Station 
Correction 

ANMO 

10610 

14 

-0.01 

ANTO 

3742 

10 

0.00 

BCAO 

7545 

6 

-0.02 

CHTO 

3888 

10 

0.08 

GRFO 

4706 

13 

-0.16 

GUMO 

7144 

15 

-0.04 

KAAO 

1888 

9 

0.10 

KONO 

4377 

15 

0.10 

MAIO 

2173 

6 

0.11 

MAJO 

4905 

14 

-0.06 

SHIO 

2926 

4 

0.07 

ZOBO 

15250 

4 

-0.09 










Figure 21. Estimated moments for the largest NTS and East Kazakh explosions since April 1976. 


Event 


Log Moment (Newton-Meters) 



































■mwspis 


moments for the largest explosions at NTS and East Kazakh since this 
date. Using a Mueller-Murphy explosion model (Mueller and Murphy, 
1971), the log moment of a 150 kiloton explosion in a tuff/rhyolite 
medium should be 16.24. This agrees very well with the upper limit of 
the recent NTS explosion estimated moments. 

A surprising result is the low values of moment estimated for 
the East Kazakh explosions. None of the explosions has an estimated 
log moment greater than 16.0. The upper limit of the log moments for 
East Kazakh is about 0.2 less than the upper limit of log moments for 
NTS. This contiasts dramatically with the NEIS m^ values for the 
explosions (also shown in Figure 21). Several of the East Kazakh 
explosions have m^ greater than 6.1, while none of the NTS explosions 
has an m^ greater than 5.8. Comparison with theoretical explosion 
models makes the contrast even larger. It is commonly believed that 
the source medium at the East Kazakh test site is a hard "granite- 
like 11 material (e.g., Alewine and Bache, 1983), but the log moment for 
a 150 kiloton Mueller-Murphy explosion in granite is 16.52, more than 
0.5 greater than the largest estimated moment at East Kazakh. 
Furthermore, the difference between m^, and M s due to local elastic 
properties at the source should be larger for a low velocity material 
then for a high velocity material (Stevens and Day, 1985). 

Several factors, including tectonic strain release, differences 
in body wave attenuation, or an unusual source medium could account 
for these differences. Correction for tectonic strain release could 
increase the moments of some events to the level of the largest NTS 
explosions (Stevens, 1982; Sykes and Cifuentes, 1984), but could not 
increase them to the level expected for a 150 kiloton granite source. 
The East Kazakh explosions are therefore not consistent with a 
"granite-1 ike” source medium. One possible explanation for the large 
body waves and small surface waves from East Kazakh explosions would 
be the presence of a source medium at the test site that resulted in a 
highly peaked source spectrum with a large amplitude near one Hz and 
an unusually small long period level. 
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X. THE EFFECT OF SOURCE REGION STRUCTURE ON MOMENT ESTIMATES 


It is important to assess the consequences of the constraints 
and assumptions that were used in making our moment estimates. The 
two most important constraints were the model for the elastic 
structure of each source region and the constraint placed on the upper 
mantle Q. In this section and the following section, we examine the 
variations in moment estimates that result from changes in these 
assumptions. 

The East Kazakh source region structure is based on the best 
available information about the region, but this information cannot be 
regarded as definitive so we want to estimate the maximum reasonable 
error in our excitation function for this region. In Figure 22, we 
show the excitation functions S} calculated from all of the earth 
structures that we have found by surface wave inversion. There is a 
substantial amount of variation in the excitation functions with a 
total variation of about a factor of two. The smallest excitation 
functions, however, belong to oceanic and mixed paths, and the largest 
are for thick crust, low velocity paths to the southeast of East 
Kazakh. In Figure 23 we show the excitation functions for all of the 
continental paths out of East Kazakh. The variation in is reduced 
substantially. The lowest amplitude excitation function belongs to 
the high velocity shield path from East Kazakh to station KONO. 

To estimate the uncertainty in the excitation function, we used 
the East Kazakh to KONO structure, with the surface velocities 
modified to be the same as in the structure used for our Q/moment 
inversions described earlier. We calculated the excitation function 
for this structure and inverted the spectra for all paths out of East 
Kazakh to obtain new moment estimates. The log moments were increased 
by 0.04 on the average. 
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Source Region Excitation Function (SI) - All Paths 
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Figure 22. Excitation function Sj calculated using the structures found for all paths 
East Kazakh. 



SI for Continental Paths from East Kazakh 
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igure 23. Excitation functions calculated using structures for all continental paths out of East 
Kazakh. 









The most important feature of the NTS source region structure 
is the shallow shear velocity. In Figure 24 we show the variation of 
the source region excitation function calculated at a frequency of 
0.05 Hz as a function of the shear velocity at the source depth of one 
kilometer. In these calculations, the shear velocity in the upper 1.5 
kilometers of the NTS source region model was varied, the 
compressionaI velocity was set using a constant Poisson’s ratio of 
0.27 and the density in this layer was held fixed. The phase and 
group velocities at this low frequency depend only weakly on the shear 
velocity, but the.excitation function S^ is quite sensitive to the 
shallow velocity structure, varying approximately as Sj ~ /T®^. 

The shear velocity at the source depth used in the earlier 
moment estimates for NTS was 2.0 km/sec. Bache, et aL. (1975) give 
shear velocities in tuff and rhyolite at NTS ranging from 1.5 to 2.5 
km/sec, so this is a good average value for the shear velocity at NTS. 
The lowest hard rock velocities at NTS are found in dry tuff, which 
have shear velocities of about 1.0 km/sec. To put an upper bound on 
the NTS source region excitation function, we modified the NTS model 
to have a shear velocity in the upper 1.5 km of 1.0 km/sec, calculated 
the excitation function for this structure, and reinverted all of the 
spectra for NTS paths to obtain new moment estimates. Log moment 
estimates were reduced by 0.08 on the average. 

The excitation functions for NTS and East Kazakh are shown in 
Figure 25 together with the excitation functions calculated from the 
two perturbed structures described here. The excitation functions are 
very frequency dependent, so it is not possible to tel I at a glance 
what the effect of the excitation function is on the moment estimates. 
To estimate the effect of the different source regions used for NTS 
and East Kazakh, we re inverted the spectra for East Kazakh paths using 
the NTS source region model. The log moments were decreased on the 
average by 0.02. 
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20 Second Sensitivity to Shollow Shear Velocity 



Shear Velocity (m/sec) 















Low Velocity NTS 



ixcitation functions calculated from NTS and East Kazakh source region structures. The 
lower dashed line is the excitation function derived from the East Kazakh to KONO path 










Correction for source region structure reduces estimated 
moments for NTS explosions by approximately 0.02 relative to estimated 
moments for East Kazakh explosions. This correction would be larger 
if either the NTS near surface velocities were lower than the values 
used, or if the East Kazakh crustal velocities were higher than those 
used in the model. Note also that the correction to the spectral 
amplitudes increases with frequency. The correction to a 20 second 
spectral magnitude should be approximately 0.09. The correction to 
moment is smaller because the moment estimate depends more strongly on 
long period amplitudes than on short period amplitudes, and because 
there is a tradeoff between the source region correction and the 
attenuation correction. For a given observed Rayleigh wave spectrum, 
moment estimation using the NTS structure will result in a larger 
attenuation correction than moment estimation using the East Kazakh 
structure, since the inversion will assume (properly) that if the 
observed spectrum had contained the amount of high frequency energy 
shown in the NTS excitation function in Figure 25 at the source, it 
must have been attenuated by the time it reached the receiver. 

One final point should be made about the source region 
correction. In Section II, we defined Sj so that we find the 
normalized moment Mg rather than the moment Mg- This choice is 
appropriate since Mg is a representation of the source of the surface 
waves that is independent of local material properties, and since 
there is no evidence that Mg correlates with yield better than Mg. We 
can, of course, convert Mg to Mg if local velocities are known through 
the relation 

M ,sL M ' 

0 3„ 20 

Since Poisson’s ratio is probably lower at East Kazakh than at NTS 
(because of stronger, more competent rock), this conversion would 
increase estimated moments at NTS relative to estimated moments at 
East Kazakh. 
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XI. THE EFFECT OF MANTLE Q CONSTRAINTS ON MOMENT ESTIMATES 


A good check on the validity of the mantle Q constraint is the 
variation of the moment estimates as a function of distance. In 
Figure 26, we show the station residuals for all NTS and East Kazakh 
paths plotted versus distance. The station residuals are small at 
close range, approximately 0.1 for most stations less than 5000 
kilometers from a test site. The magnitude of the residuals increases 
with distance, as is to be expected with the increasing complexity of 
the path at large distances, but the residuals scatter evenly around 
zero. The lack of any trend in the data suggests that the attenuation 
correction has been properly estimated with the constraint of /J/Q = 40 
between 120 and 150 kilometers depth. 

In Figures 27, 28, and 29 we show station residuals with 
constraints of /J/Q =25, 50 and 100 respectively. With the lower 
attenuation of /?/Q =25 (Q ~ 180), the station residuals clearly 
increase with distance indicating that the moments are undercorrected 
for attenuation. The relative moments at the two test sites are not 
changed significantly, however. The log moments are decreased by 0.05 
at NTS and by 0.07 at East Kazakh. Increasing attenuation slightly 
with P/H - SO (Q «* 90), the station residuals show little trend with 
distance and the moments are increased by 0.03 at East Kazakh and by 
0.04 at NTS. If attenuation is increased to pi Q = 100 (Q ~ 45), then 
the station residuals show a strongly decreasing trend with distance 
indicating that the moments have been overcorrected for attenuation. 
In this case the moments are increased by 0.22 at NTS and by 0.23 at 
East Kazakh. 

Based on these tests, we conclude that a constraint of 
/?/Q s 40-50 in the upper mantle is appropriate. The relative moments 
between the two test sites is insensitive to the choice of constraint 
as long as the same constraint is used for both test sites and the 
distribution of stations is comparable for both test sites. 
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Figure 26. Station residuals with a mantle Q constraint of /J/Q = 40. 









Station Residual Beta/Q = 25 



Figure 27. Station residuals with mantle Q constraint of /J/Q = 25. 















28. Station residuals with a mantle Q constraint of /?/Q 











Station Residual Beta/Q = 100 



Figure 29. Station residuals with a mantle Q constraint of /l/Q 











The value of /J/Q in the mantle is not, of course, constant 
throughout the earth. Because of the finite frequency range of the 
data, we cannot directly determine the mantle Q, and our approach has 
been to use an average value as a constraint and to accept some 
scatter in the results, rather than to try to make uncertain 
corrections for variations in mantle Q. This is particularly 
appropriate for long paths, where the Q co’rrection is most important, 
since many different regions of the earth are sampled by the Rayleigh 
wave between the source and receiver. 

An interesting question is the effect of a small region of the 
travel path with a mantle Q significantly different from the average 
value. As we mentioned earlier, for example, interstation Q studies 
have found much lower mantle Q values near NTS (e.g., Lee and Solomon, 
1979) than near East Kazakh (e.g., Patton, 1980). If the mantle Q is 
lower than our average value, then the attenuation coefficients will 
be underestimated. From Equation ( 6 ), we find for a moment estimate 
M 0 made from an observed spectrum u(«), 

AlnM 0 = A(r 7 («)) (9) 

so if 7 is underestimated, Mg will also be underestimated by an amount 
that depends on the length of the travel path with the different Q 
value. Comparing results from our test cases above, we find that 
changing mantle /J/Q from 40 to 25 decreases 7 by approximately 
0.3 x 10~ 4 /km over the frequency range 0.02 to 0.05 Hz. Changing pf Q 
from 40 to 100 increases 7 by approximately 1.0 x 10" 4 /km over the 
same frequency range. If p /Q = 100 were a better value for mantle I) 
for a region approximately 1000 km around NTS, then our log Mg 
estimates for NTS explosions would be low by about 0.04. If /7/Q = 25 
were a better value for mantle Q for a region 1000 km around East 
Kazakh, then our log Mg estimates would be high by approximately 0.01. 
Making this additional correction would therefore increase the NTS 
moments relative to the East Kazakh moments by approximately 0.05. 
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SUMMARY AND CONCLUSIONS 


In summary, the most important conclusions of this report are: 

* The path correction procedure does provide stable and 
accurate measurements of surface waves. 

* Typical network scatter in log moment is about 0.1 
compared to 0.2 or more for network scatter in M $ . 

* Dispersion is not a problem since this is a spectral 

method. Data from stations at close range (less than 

30 degrees) can be used without difficulty. 

* The spectrum is flattened by application of the path 

correction, so the results are nearly frequency 
independent. Moments are estimated from an average 
over a range of frequencies. 

* Little variability in spectral amplitudes results from the 
elastic structure of the travel path. 

* The standard deviation on the logarithm of the path 

amplitude function over all of the paths processed is 
only 0.025. The extreme values differ by only 0.15. 

* This result contrasts dramatically with time domain 
amplitudes which have large variations due to path 
dependent dispersion. 

* The elastic structure of the source region can 
significantly affect surface wave amplitudes. 

* Assuming that the surface velocities are higher at East 
Kazakh than at NTS, surface wave amplitudes should be 
larger for NTS explosions than for East Kazakh 
explosions with the same moment. 
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* Correction for source region structure based on the 
best models available for the two test sites reduces 
the estimated log moment of NTS explosions by 0.02 
relative to log moments for explosions at the East 
Kazakh test site. 

Attenuation differences can cause large variations in 
surface wave amplitudes at distant stations. An important 
check on the attenuation correction is the variation of 
moment estimates with distance. 

* When the Q inversion is performed with the mantle Q 
constrained to a value of approximately 100 at depths 
of 120 to 150 kilometers, the moment estimates remain 
constant (on the average) with distance. 

* When the Q inversion is performed with a Q constraint a 
factor of two higher or lower, the moment estimates 
show a significant trend with distance indicating under 
or overcorrection for attenuation. 

* The relative moments between the two test sites are 
insensitive to changes in the Q constraint. 
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